JOURNAL OF GEOPHYSICAL RESEARCH, VOL. 102, NO. A4, PAGES 7509-7521, APRIL 1, 1997 


NASA-CR-204685 

Photoelectron effects on the self-consistent potential in the 
collisionless polar wind _ r 

G. V. Khazanov, M. W. Liemohn, and T. E. Moore - ^ 

Space Sciences Laboratory, NASA Marshall Space Flight Center, Huntsville, Alabama 

Abstract. The presence of unthermalized photoelectrons in the sunlit polar cap leads to an 
enhanced ambipolar potential drop and enhanced upward ion acceleration. Observations in the 
topside ionosphere have led to the conclusion that large-scale electrostatic potential drops exist 
above the spacecraft along polar magnetic field lines connected to regions of photoelectron 
production. A kinetic approach is used for the 0 + , H + , and photoelectron (p ) distributions, while a 
fluid approach is used to describe the thermal electrons ( e ) and self-consistent electric field (£| ( ). 

Thermal electrons are allowed to carry a flux that compensates for photoelectron escape, a critical 
assumption. Collisional processes are excluded, leading to easier escape of polar wind particles and 
therefore to the formation of the largest potential drop consistent with this general approach. We 
compute the steady state electric field enhancement and net potential drop expected in the polar 
wind due to the presence of photoelectrons as a function of the fractional photoelectron content and 
the thermal plasma characteristics. For a set of low-altitude boundary conditions typical of the 
polar wind ionosphere, including 0.1% photoelectron content, we found a potential drop from 500 
km to 5 R e of 6.5 V and a maximum thermal electron temperature of 8800 K. The reasonable 
agreement of our results with the observed polar wind suggests that the assumptions of this 
approach are valid. 

lighter ions must be dragged away from the upper atmosphere 
by the escaping photoelectron flux. Lemaire [1972] pointed 
out that the additional electrostatic force acting on the ions 
because of the escaping photoelectrons can accelerate the ions 
to higher velocities as long as the photoelectron flux is larger 
than the escape flux of the thermal electrons. This idea was 
developed further by Barakat and Schunk [1984] by including 
precipitating hot electrons of magnetospheric origin in the 
polar cap region. Several tutorials by Schunk [1986, 
1988a,b] and the recent review by Ganguli [1996] provide a 
complete picture of the historical development of polar wind 
studies. 

Recently, motivated by Akebono satellite measurements, 
Tam et al. [1995a] have developed a numerical steady state 
polar wind solution that is continuous from the subsonic col- 
lisional regime at 500 km to the supersonic collisionless 
regime at 2 Rg. Wilson et al. [1996] reexamined this problem 
in order to reduce the large acceleration of 0 + ions to super- 
sonic speeds at altitudes below 600 km and to reduce the high 
thermal electron temperatures in excess of 40,000 K at 1000 
km that came from the results of the Tam et al. [1995a] model. 
They found that by balancing the photoelectron flux with a 
downward electron flux, instead of an upward 0 + flux, a poten- 
tial drop of only 5 to 6 eV results between 500 km and 3 R E . 
They also introduced a precipitating magnetosheath electron 
flux, or polar rain, to help balance the photoelectron outflow. 

The escape of photoelectrons from the sunlit polar cap 
causes a potential drop to develop that leads to additional ion 
acceleration. Observations in the topside ionosphere by ISIS 
1 [Winningham and Heikkila , 1974] and DE 2 [ Winningham 
and Gurgiolo , 1982] have led these authors to conclude that 
large-scale electrostatic potential drops exist above the space- 
craft along polar magnetic field lines connected to regions of 
photoelectron production. Pollock et al. [1991] examined the 
occurrence rate and magnitude of Field-aligned electrostatic po- 


1. Introduction 

The polar wind, ion upwelling events, and plasmaspheric 
refilling are major sources of mass, momentum, and energy for 
the entire magnetosphere. As the plasma flows up and out of 
the topside ionosphere, the flow conditions change from sub- 
sonic to supersonic, from collision-dominated to collision- 
less, and from 0 + dominance to H + dominance. In the colli- 
sionless regime, the ion velocity distributions become highly 
non-Maxwellian, and the coupling between various plasma 
species occurs through the development of a self-consistent 
potential. 

The reason for the formation of a self-consistent potential 
in the collisionless plasma is quite clear. High-mobility elec- 
trons tend to overtake ions. As a result, the electric neutrality 
of the plasma is violated and an electric field appears which 
constrains the electrons, forcing them, on average, to travel 
together with the ions. This field also significantly affects 
the motion of the ions by accelerating them. The electric 
force eE is proportional to the electron temperature T e and has 
to be combined with the gravitational force mg to determine 
the pressure distribution of the ions. 

Photoelectrons, which form due to ionization of the atmo- 
sphere by solar radiation, can alter the self-consistent poten- 
tial in the space plasma. The presence of the enhanced high- 
velocity tail in the electron distribution will increase the 
number of fast ions. Because of enhanced ion acceleration in 
an expanding plasma, an initial superthermal electron distri- 
bution function may be changed. 

The polar wind concept was introduced by Banks and Holier 
[1968] and Axford [1968]. The latter also suggested that the 

This paper is not subject to U. S. copyright. Published in 1997 by the 
American Geophysical Union. 

Paper number 96JA03343. 


7509 



7510 


KHAZANOV ET AL.: PHOTOELECTRON EFFECTS IN COLLISIONLESS POLAR WIND 


tential drops over the ionospheric polar caps. In that study, 
signatures in upgoing and downgoing photoelectrons were 
measured in the topside ionosphere using the low-altitude 
plasma instrument (LAPI) on Dynamics Explorer 2 (DE 2) 
[Winningham and Gurgiolo, 1982]. These data are compared 
with ion data obtained at high altitudes using the retarding ion 
mass spectrometer (RIMS) on the DE 1 spacecraft. Data were 
selected from intervals when DE 1 and DE 2 were approxi- 
mately along the same polar cap magnetic field line and when 
upflowing 0 + beams were observed in the RIMS data. Pollock 
et al. [1991] presented one case in which the comparison of 
data from the two DE spacecraft is quite favorable regarding its 
interpretations in terms of a field-aligned electrostatic poten- 
tial drop. 

What is the contribution of photoelectrons in the formation 
of a field-aligned electrostatic potential drop over the polar 
cap? What is the largest potential drop that should be ex- 
pected due to the presence of photoelectrons in the polar wind? 
What is the quantitative relation between the thermal plasma 
and photoelectron parameters that leads to certain values of 
the self-consistent electric field in the polar wind? How do the 
thermal plasma parameters affect the result for a given photo- 
electron density? These questions are the primary motivation 
for the study presented in this paper. 


2. General Relations 


In this paper, we will consider the calculation of the upper 
limit of the potential drop that should be expected in the polar 
wind due to the presence of photoelectrons. The exclusion of 
collisional processes leads to easier escape of polar wind par* 
tides and therefore to the formation of the largest possible 
field-aligned electrostatic potential drop over the polar cap. 
This argument is consistent with the discussion of Tam et al. 
[1995b] and should only be valid above the region of the 
source where transport is the dominant process in polar wind 
formation. In this paper, a kinetic approach is used for the 0 + , 
H + , and photoelectron ip ) distributions, and a fluid approach is 
used to describe the thermal electrons (e) and self-consistent 
electric field (£ g ). 

The steady state collisionless kinetic equation can be pre- 
sented in the following form: 
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where v is velocity, p is the cosine of pitch angle, 
fa - faf^fipS) * s the distribution function of species a, s is 
distance along the geomagnetic field line, B is the geomag- 
netic field, and g is the gravitational acceleration. 

The most delicate part of the polar wind problem in the 
presence of photoelectrons is photoelectron-thermal electron 
interactions. In the absence of collisions, the electron plasma 
component is coupled with 0 + and H + ions based on the quasi- 
neutrality and currentless conditions through the self-consis- 
tent electric field. In a steady state collisionless plasma, mass 
and energy conservation equations for the total electron com- 
ponent can be presented in the form 


d_ 

ds 
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where n t w, and Q are the moments of electron distribution 
function 


V 

«(j) = — f/!v/(v,/i,j)rf 3 v (4) 

G( j )=yj/ iv3 /( v .^. i y 3v 

v 

Here d*v = 2jcw 2 d\dp is the volume in the velocity space. 
The term with gravitation in (3) can be omitted for electrons. 

The electron population in (4) can be separated into two 
parts: thermal electrons, with distribution function f e \ and 
photoelectrons, with distribution function f p . Accordingly, 
the mass and energy conservation equations can be presented 
as 



Using the kinetic equation (1) for the photoelectron distribu- 
tion function (a=p), we can find mass and energy conservation 
for this portion and exclude the photoelectron component 
from (5) and (6). This leads to the set of thermal electron mass 
and energy equations 



As we pointed out before, these moment equations should be 
coupled with 0 + , H + , and photoelectrons based on the quasi- 
neutrality and currentless conditions 

*e M + n p (■*) = *o + ( S ) + n H+ (*0 


ne(s)u e (s)+n p (s)u p (s) = n Q <r (s)u Qt (f) + n H . (j)u h * M ( 9 ) 

where densities and fluxes in (9) for 0 + , H + , and photoelec- 
trons must be found based on the solution of the kinetic 
equation (1). 


3. Solution of the Kinetic Equation 

The solution of the steady state collisionless kinetic 
equation for polar wind applications has been established and 
discussed in great detail in several papers by Lemaire and 
Scherer [1970, 1971, 1972]. They classify the polar wind par- 
ticles into four categories [Lemaire and Scherer , 1972]: (1) 
ballistic particles that cannot escape; (2) particles with 
enough energy to escape to interplanetary space; (3) trapped 
particles with two mirror points along the field line; and (4) 
incoming particles from the outermost regions. 
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For the purposes of our study, it is assumed that the incom- 
ing particles are absent. This assumption leads to an electro- 
static potential distribution with no possibility of trapping 
particles, and so this population is also omitted (see section 7 
and Figure 1). 

Since the plasma is considered to be collisionless, the time- 
independent distribution function in (1) depends only on the 
particle's total energy 


reach s' and contribute to the particle density rt a (.y'). 
Particles having 

E<e a U{s’)-m a H'{s’) + M a B{s') (16) 

are reflected before reaching s'. Using these as limits of inte- 
gration in (13), the particle density and flux for species a can 
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be written as 


(io y— 


„ \ 3/2 °o 

-ftr) Jr 


and magnetic moment 
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The potential differences U(s) and •f'U) can be related to the 
parallel electric and gravitational fields from the following 
relations 


x J dE 
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{E-e a V+m a 'i'-MB) U1 


n a u a = 
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^ 2 ) and depend on the electric and gravitational potentials, mag- 
netic field, and ionospheric distribution function. 


The particle distribution is constant along any particle tra- 
jectory which is characterized by the constants of motion; 
therefore the distribution function at any point along the ge- 
omagnetic field can be calculated provided it is known at some 
boundary in the ionosphere above the region of the source. 
For convenience, we will transform our distribution function 
and integration variables from (v, p) to the constants of 
motion (E, M). The particle density and flux for species a can 
then be written as 


\Uejm 

2 V w cr J (E-e a U + m cc l F- 


MB) 


1/2 


kB 

n a u a = — 



(M.E)dEdM 


(13) 


which is similar to that found by Whipple [1977] and used by 
Miller and Khazanov [1993] to calculate a steady state solu- 
tion for the self-consistent electrostatic potential due to 
trapped plasma in the magnetosphere. 

The particle density at arbitrary positions along the mag- 
netic field can now be easily determined provided the distribu- 
tion function is known at the lower, ionospheric boundary. 
The region of determination for the integrals (limits of inte- 
gration) in (13) is obtained by considering the condition of 
reflection, which occurs for V|j =0, and is given by 


( o V /2 

V|l =f-2-J {E-e a U + m a 'F-MB)' /2 = 0 (14) 


which is a function of the total energy of the particle, the local 
electric and gravitational potentials and magnetic field, and 
the sign and the mass of the charged particle. Equation (14) 
depends only on the final position s' and is a necessary but not 
sufficient condition to guarantee accessibility of an iono- 
spheric particle to s'. For particles to gain access to s', they 
must not be reflected at intermediate positions (0 < s < s') 
along the magnetic field. According to (14), particles with 
energy 


4. Moments of the Velocity Distribution Function 

In order to calculate the hydrodynamic moments (4), the dis- 
tribution function at the lower boundary of the simulation do- 
main must be specified. For the sake of comparing our results 
with Tam et ai [1995a], we start at 500 km (low-altitude 
boundary) with upper-half Maxwellian ion and photoelectron 
distributions. It must be noticed, however, that analytical and 
numerical analyses for the polar wind [Khoyloo et al . , 1991] 
indicate that the adoption of a Maxwellian for an expanding 
collisionless plasma in the diverging geomagnetic field is in- 
consistent and causes a discontinuity in the kinetic solutions. 
Additionally, choosing the lower boundary above the plasma 
source region introduces additional uncertainties into the 
problem and predetermines the solution above the lower 
boundary level (T. I. Gombosi, private communication, 
1992). 

In our calculations, we separated particles that started from 
the lower boundary and reached point s on along the field line 
(part of the escaping population discussed by Lemaire and 
Scherer [1972]) from particles that have been reflected below 
this point (a component of the ballistic particles). In accor- 
dance with this, we call them the transient and reflected popu- 
lations. Depending on the electrostatic e a U(s) and gravita- 
tional ma^is) energies of polar wind particles, two cases 
should be considered separately; case 1, e a U(s)> m a i F(s) 1 
and case 2, e a U{s)<m a ^(s). 

Case 1 only has a transient population (denoted tr) and the 
following expressions for the density and particle fluxes 


n'a( s ) = n Oae Za 



/ 

yexp 
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JiZa) 
Bo y 2 ) 



n a ( s ) u a ( s ) = ^Oa^Oa ( |g ) 

Here subscript ”0" refers to the lower boundary level; <P is the 
error function, 


<J>(x) = (2 /Vz)f Q e ,2 dt 


E > e a U{s')~ m a '{'{s') + M a B{s') 


(15) 
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and the z and y parameters are defined as 


Za — 
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where T a is the characteristic temperature (energy) of the 
Maxwellian for species a at the base. The total density and 
particle flux that should be used in the quasi-neutrality and cur- 
rentless conditions (9) in this case are 


"a( s ) = »i'W 


(19) 


Case 2 includes both the transient and reflected (denoted ref) 
particles, and the following expressions for the density and 
particle fluxes can be found to be 
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where F is Dawson's integral, 

. . i ex 2 

F(x) = * _Jf J ** 

As was pointed out above, the reflected population in (21) 
corresponds to the particles that have been reflected below 
point s back to the base. These quantities are needed at point 
$, though, and the changing cross-sectional area of the flux 
tube must be taken into account. Therefore the total density 
and particle flux that should be used in the quasi-neutrality and 
currentless conditions (9) are 

«a (*)»<*(*) = «£(*)“£(*) (22) 

-tkVKVK'ow'w] 

Here we add to the transient particles that reached point s the 
component of particles reflected between point 5 and the upper 
level of the simulation domain s kb , The flux equation in (22) 
can be rewritten in a more simple form as 
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where y & is taken at the upper boundary level. 


5. Thermal Electron Fluid Equations 

In the previous section, the analytical expressions for the 
ion and photoelectron moments in the self-consistent 
coupling with the thermal electron fluid equations (7) and (8) 
are found as functions of the magnetic and gravitational fields 
and the self-consistent electrostatic potential. Now we should 
close the loop and find the corresponding expression for the 
density of the thermal electrons. The implicit relation for the 
density of the thermal electrons as a function of the self-con- 
sistent electrostatic field can be found from the thermal elec- 
tron momentum equation 


ds 


d[Am e n e u 2 ) 


ds 


-enE\\ = 0 


(23) 


where A is the cross-sectional area of the flux tube. Simple 
integration of this equation leads to the following expression 
for the thermal electron density 


n e (s) = 


1 


T e (s) + m e u?(s) 
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ds' > 


Now, using the condition of quasi -neutrality in (9), an itera- 
tive procedure can be used to find the distribution of the self- 
consistent potential in the polar wind plasma. Only one pa- 
rameter in (24) still remains unknown, the electron tempera- 
ture T e (s). To find it, let us present (8) in the conventional 
form for the electron temperature by separating the total ther- 
mal electron heat flux Q into the thermal energy flux Qj and 
the convective particle motion. Then, using (6) and consider- 
ing the thermal energy flux as a thermal conductivity flux 
according to Banks and Kockarts [1973], 


n T _ v t 5/2 dT e 

Q'-x° T ‘ 


(25) 


the equation for the thermal electron temperature can be writ- 
ten in the form 





XoT? 2 3T « 
B ds 


= 0 


(26) 


Here U is potential difference at the upper level of the sim- 
ulation domain, and Xo is the appropriate constant associated 
with the thermal conductivity. To calculate the thermal elec- 
tron particle flux in (26), n e u ey the currentless condition in 
(9) should be used. 

Equation (26) can be integrated twice with respect to s and 
presented as 


7*7/2 _ r 7/2 , 7 f j n e u e 
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The square brackets in (27) are taken at the upper boundary of 
the simulation domain. 

It should be pointed out here that in order to provide an ana- 
lytical description of the electron temperature, we greatly 
simplified the calculation of the electron thermal flux (25) by 
using Spitzer's expression for the conductivity coefficient. 
We will address this issue later in section 7. 

Thus our polar wind model uses the kinetic approach to cal- 
culate the distributions of photoelectrons and the 0 + and H + 
ions in the polar wind plasma as an explicit analytical func- 
tion of the self-consistent electrostatic potential (see section 
4) and an implicit analytical solution of the fluid equations 
(24) and (27) to calculate the distribution of the thermal elec- 
trons and their temperature, which also contains the self-con- 
sistent electrostatic potential. Because of the quasi-neutrality 
and currentless conditions (9), all polar wind plasma parame- 
ters and components are coupled and can be determined based 
on an iterative scheme. 


6. Iterative Solution Method 

With an initial guess for the U(s) and T e (s) distributions, 
it is possible to iterate to a convergent solution that satisfies 
all of the equations. At each iteration step, the quasi-neutral- 
ity condition is satisfied using a root-finding method at a 
given altitude point, and then the electron temperature is cal- 
culated based on this new potential and electron density. We 
start each iteration at the base altitude and work up the field 
line so the integrals for n e and T e are consistent with the cur- 
rent iteration. The upper boundary terms are taken from the 
previous iteration, and the process continues until U(s) and 
T e (s) converge at every step in the spatial domain. 

The quasi-neutrality condition in (9) can be rewritten as 

F[V) m '2 l *a*a(U) (28) 

a 

and then solved for F(U)=0. Brent’s method was used for this 
study, which combines the speed of inverse quadratic interpo- 
lation with the reliability of the bisection method, and is 
guaranteed to find the root given an interval that contains one 
[Press et al. , 1992]. After finding an interval of U that con- 
tains a root of F, this method finds the exact U that balances 
all of the densities, no matter how rapidly the densities change 
with U. 

Since the moments of the kinetic equation only depend on 
the current spatial step and the upper and lower boundaries, 
any altitude grid could be used for these equations. However, 
the numerical integrations in n e and T e dictate that a reason- 
able step size must be used. This is discussed further in section 
7.2. 

One possible problem can arise from the relation between 
the electrostatic and gravitational potential differences, eU 
and mY. The kinetic formulae above are derived assuming that 
each plasma species always satisfies either case 1 or case 2. 
This is not a problem for H + and the photoelectrons, where the 
particle mass is sufficiently small so eU always dominates. 
For 0 + , however, the electrostatic and gravitational potential 
differences can be similar. For instance, ^(3 Re) is 7 eV, 
while eU(3R E ) has been shown to range from 6 eV [Wilson et 
al. t 1996] to 12 eV [Tam et al ., 1995a]. Thus it is unknown 
what case 0 + will follow, and the potentials could even cross 
somewhere inside the spatial domain. In general, if m were 
to dominate at low altitudes and then be surpassed by eU 


somewhere up the field line, the governing equations would 
switch from case 2 to case 1 . In the case 2 equations, the upper 
boundary would be this crossover altitude, and similarly, in 
the case 1 equations, this spatial point would be the lower 
boundary. Should this type of crossover occur in our calcula- 
tions, the crossover boundary is determined from the previous 
iteration, but this point is corrected as the calculation pro- 
ceeds up the field line, if necessary. This could make the lower 
altitude points inconsistent with the higher altitude solutions, 
but since the distribution of U is converging with each itera- 
tion, these inconsistencies are temporary and will not affect 
the coherence of the final result. In summary, our algorithm 
checks for the condition of crossing potentials in the 0 + cal- 
culation and handles it accordingly. 

7. Results 

As we mentioned in the Introduction, the purpose of this 
study is to answer the following questions: What is the con- 
tribution of photoelectrons in the formation of a field-aligned 
electrostatic potential drop over the polar cap? What is the 
largest potential drop that should be expected due to the pres- 
ence of photoelectrons in the polar wind? What is the quanti- 
tative relation between the thermal plasma and photoelectron 
parameters that leads to certain values of the self-consistent 
electric field in the polar wind? How do the thermal plasma pa- 
rameters affect the result for a given photoelectron density? 
These will now be addressed by iteratively solving equations 
(9), (27), and (28) along a polar cap field line from 500 km to 
5 R e . 

Although there are other factors that affect polar ion mo- 
tion, such as currents, precipitation, and plasma waves, our 
objective is to show the impact of just one factor, photoelec- 
trons. Since we are solving the collisionless kinetic equation, 
we believe we should obtain an upper limit to the photoelec- 
tron influence for a given set of input parameters. 

The initial guess for U(s uh ) and T e (s uh ) used depends on n p0> 
the photoelectron content at the base, and ranged from U(s ub ) 
=3 V and T e (s ub ) =2000 K for n p(? =0.0% up to U(s ub ) =7 V and 
T e (s ub ) =25000 K for n p0 = 1.0%. Equations (24) and (27) are 
sensitive to the U and T e distributions, and so these initial 
guesses must be reasonable for the solution to converge to a 
solution. If U/T e is too large, the initial n e distribution will be 
very low, resulting in large, unrealistic velocities that can 
cause (27) to have a negative integrand. It was found, how- 
ever, that the result is independent of the initial guess at 
U(s ub ) and T e (s ub ), as long as they produce a positive inte- 
grand in (27) for the first iteration. 

Several initial spatial distributions of these two parameters 
were also tested. Such choices for U(s) as a linear fit, U(s) 
proportional to [5/(5 M ^-5o)] 1/2 , and U(s) proportional to s ) 

produce distributions that converge to the same solution. 
Likewise for T e (s), a linear fit, T e (s) proportional to [$/(j w *- 
j 0 )] l/2 , and TJs) proportional to sin[.5ns/(s ub -s 0 )] all produce 
distributions that converge to the same solution. Thus the ini- 
tial spatial distribution of U and T e is irrelevant, as long they 
produce a positive integrand in (27) for the first iteration. 

As follows from the equation for electron temperature (27), 
we need two boundary conditions to calculate T e \ a lower 
boundary temperature T € q and a temperature gradient at the top 
of the simulation domain [dTg/ds] s ub. As we mentioned before, 
the purpose of this study is to investigate the role of photo- 
electron escape on the thermal plasma parameters. Because 
photoelectrons form in the ionosphere and flow out along the 
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geomagnetic field lines, their influence on the thermal plasma 
in the collisionless limit is only through the formation of the 
self-consistent potential, as seen in (27). To focus on this ef- 
fect, we will omit magnetospheric energy fluxes by assuming 
[dT /ds] s ub=0. We will show a case, however, with a source 
term in the electron temperature equation that accounts for 
Coulomb collision heating of the thermal electrons by the 
photoelectrons (see section 7.3). 

7.1. Photoelectron Density Dependence 

To illustrate the photoelectron dependence of our results, 
typical ionospheric parameters for the polar region were cho- 
sen for the low-altitude boundary conditions. The characteris- 
tic temperature for each particle species is 2000 K for the ions, 
2500 K for the thermal electrons, and the photoelectrons have 
a characteristic energy of 20 eV. The 0 + density at 500 km is 
6xl0 4 cm' 3 , and H + has lxlO 3 cm' 3 , for a total ion density of 
n tot =6.1xl0 4 cm' 3 . The photoelectron density is varied from 

0. 00% to 1.00% of this number, and the thermal electron den- 
sity at the base is determined by quasi-neutrality. The density 
scale is n p0 = 0.01% to 1.00% on the Plates because of the log- 
arithmic axis. A photoelectron concentration of 1.00% at 500 
km is huge and unrealistic but is shown to illustrate an extreme 
upper limit to the photoelectron influence. Ionospheric polar 
cap photoelectron densities rarely exceed 0.10% [Khazanov 
and Liemohn, 1995]. 

Although the format of the Plates is unconventional, it is 
convenient to show the dependence of the results on the pho- 
toelectron concentration throughout the polar altitude range. 

The electrostatic potential (p as a function of the photoelec- 
tron relative density at the base and altitude is shown in Plate 

1. The potential difference, U y is related to the potential by 
U(s)=<p 0 -<p(s) y where <p 0 is the potential at 500 km and s ub is a 
reference altitude where the potential equals zero. For the 
boundary conditions listed above, (p 0 is 2.87 V for 
n p o= 0.01%, while the total potential for 1.00% is 7.25 V. It 
can be seen that the distribution of <p is almost constant for 
fl p o above 0.10%, with (p 0 changing less than a volt. The po- 
tential also starts to reach an asymptote at low concentrations 
*lso. For any n p o , however, most of the potential drop occurs 
within the first 2 to 3 /?£, smoothing out to zero at the top. 

Such electrostatic potential behavior, combined with the 
gravitational potential, will not create a trap with two mirror 
points in the spatial domain. Figure 1 shows the total poten- 
tial energy of the three species that were treated kinetically for 
a photoelectron concentration of 0.10% at the base. As you 
can see, there is no possibility of a "potential well" to trap 
these particles along the field line. This justifies our assump- 
tion earlier that we do not need to include a trapped population 
in our model. 

Plate 2 shows the thermal electron temperature for the same 
calculation. While T € q is 2500 K for all photoelectron densi- 
ties, the temperature at the top varies from less than 3600 K at 
0.01% up to 17,300 K at 1.00%. At any given n p0t T e rises 
quickly through the first one to two R E , then flattens into a 
nearly constant temperature for the remainder of the field line. 
This is due to the thermal conductivity redistributing the en- 
ergy, smoothing the temperature distribution from the low-al- 
titude source region into the upper altitudes. The energy source 
for this rise in T e is from the electrostatic coupling of the pho- 
toelectrons with the thermal plasma. Notice that the maxi- 
mum temperature for a realistic photoelectron density 


(npfl<0.10%) does not exceed 9000 K, in disagreement with 
the results of Tam et ai [1995a], but consistent with the ob- 
servations of Abe et ai [1993]. 

Another interesting result is the ratio of 0 + to H + densities, 
shown in Plate 3. The classical polar wind description has 0 + 
dominating at low altitudes and H + dominating at high alti- 
tudes. This would appear as a ratio greater than one decreasing 
to a ratio less than one, reaching unity a few thousand kilome- 
ters up the field line. This is exactly what is seen at low n p0 
concentrations, with a ratio of unity being reached at less than 
an R e above the surface and then H + completely dominating 
the total ion density. As we move up in photoelectron con- 
centration, though, the situation changes as the electrostatic 
potential grows large enough to drag the 0 + ions against grav- 
ity to higher altitudes. By n p o=0AQ% y the two ions have sim- 
ilar densities along most of the field line, and, above this pho- 
toelectron percentage, 0 + dominates over the entire spatial 
region. 

The 0 + /H + density ratio increases near 5 R E at large n p0 
concentrations. The ratio at the base is 60 and, for 1.00% 
photoelectrons, drops to 28 near 2 R E and then increases up to 
72 at 5 R E . This behavior can be explained by examining z(s) 
for 0 + , which describes the relationship between the gravita- 
tional and electrostatic potential differences. At 1.00%, z is 
zero at the base, increases to 4.1 at 2 R E> and then decreases to 
3.65 at the top. Since exp(-z) is a multiplier in the 0 + density 
calculation (see equations (20) and (21)), this increase below 2 
R e means that n 0 + decreases more rapidly than n H +, while the 
decrease in z above 2 R E indicates that no+ decrees less 
rapidly than n h+. These trends in the densities are reflected in 
the ratios shown in Plate 3. 

Another feature that should be discussed is the ion speeds, 
illustrated here with the Mach number. Mq+ is given in Plate 4 
and A/ h + is given in Plate 5. For these plots, the Mach number 
is defined as the bulk velocity at a given altitude divided by the 
thermal velocity at the base. Therefore all of these Mach 
numbers are relative to an ion thermal speed for a 2000 K 
Maxwellian distribution, that is, 0.81 km s' 1 for 0 + and 3.2 
km s 1 for H + . 

In Plate 4, Af 0 + changes from subsonic for the entire field 
line at 0.01% to supersonic for most of the field line at 1.00%. 
This subsonic region corresponds to the classical polar wind 



Figure 1. Total potential energy for 0 + , H + , and electrons as 
a function of altitude for n pO =0.10%. 
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Plate 1. Electrostatic potential as a function of photoelectron density and altitude. Low-altitude boundary 
conditions are 0 + density of 6xl0 4 cm'- 1 , H + density of IxlO 3 cm' 3 , r,o= 2000 K, T e g=2500 K, and £ p o= 20 eV. 
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Plate 2. Thermal electron temperature as a function of photoelectron density and altitude. Same boundary 
conditions as Figure 1. 
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Plate 3. 0 + /H + density ratio as a function of photoelectron density and altitude. Same boundary conditions 
as Figure 1 . 
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Plate 4. 0 + Mach number as a function of photoelectron density and altitude. Same boundary conditions as 
Figure 1. 
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Plate 5. H + Mach number as a function of photoelectron density and altitude. Same boundary conditions as 
Figure 1. 


region discussed with Plate 3, with velocities orders of magni- 
tude less than the thermal speed for the first few R E up the field 
line. In the opposite extreme, for example, n p o= 1.00%, 0 + 
becomes supersonic by 1500 km, reaches a Mach number of 
2.7 near 2 Re, and then slowly lessens to Afo+“l-l at the top. 
This appears to be an asymptotic limit for the 0 + Mach num- 
ber at high altitudes, corresponding to a velocity just under 1 
km s' 1 . This maximum in the 0 + velocity near 2 Re is again 
due to the difference in the spatial dependencies of the gravita- 
tional and electrostatic potentials. 

The protons, however, are always supersonic for the given 
initial parameters (Plate 5). Mach 5 is achieved between 2000 
and 2700 km altitude after a quick acceleration, and then it 
seems to increase more slowly above this level. Again, a clear 
separation between the classical polar wind and a different 
state is evident. At low ^values, the Mach number reaches 7 
to 8 at the top, but at higher n p o values, the speeds break Mach 
11 (35 km s' 1 ). Even with this difference at higher altitudes, 
the main acceleration region for H + is at low altitudes, regard- 
less of the photoelectron concentration. 

7.2. Boundary Condition Dependence 

It is also useful to determine the influence of changing vari- 
ous other parameters of the calculation. Increasing the number 
of spatial points has very little affect on the results. In 
doubling the number of altitude levels from 100 to 200, T e and 
<p are not significantly changed (<1% difference), showing 
that the numerical integrations are consistent. Changing the 
upper boundary level also has very little affect on the results. 
Dropping s ub to 3 R E creates less than 5% difference in T e and 
less than 1% difference in <p. 

The choice of characteristic energies at 500 km can be more 


significant, but they are still not crucial parameters. 
Increasing the photoelectron energy to 40 eV, or even 60 eV, 
causes only a minor correction to the results, increasing T e by 
roughly 15% for each of these steps. This acts to decrease n p 
and increase n t and /iq + slightly, and 0 + becomes supersonic 
earlier, dropping from 2 R E for E p0 = 20 eV to 1.5 R E for 
E p o=60 eV. The electrostatic potential is not significantly 
changed. Increasing T e0 to 4000 K only affects the first few 
spatial points, and by 1500 km, the influence on T e (s) is less 
than 1%. The variable <p 0 is 0.1 V higher, but (p(s) drops faster 
in these initial points, and is then similar throughout the spa- 
tial range. This acts to increase the 0 + velocity near the base, 
but the velocity is unchanged above 1 R E . Changing all ther- 
mal temperatures at the base to 4000 K acts to decrease (p by 
20% below 2 R E . At this altitude, hh+ ^ as doubled and riQ+ has 
tripled over the previous results. There is a similar jump in the 
velocities, but the Mach numbers are very similar to the pre- 
vious case. This decrease in the potential is because the ions 
are already moving faster, therefore less of an assist from the 
electrons is needed to lift them along the field line. 

Another point to discuss is the dependence of the results on 
the reflected particle population. Removal of these terms from 
the case 2 equations yields less than a 1% correction in the 
electrostatic potential and less than 0.1% difference in the 
electron temperature. The spatial distribution of the photo- 
electron density is affected, however, dropping to almost half 
of its value of when reflection was included. The 0 + density 
distribution is also changed, dropping by an amount matching 
the photoelectron depletion. There is no significant effect on 
H + or the thermal electrons. From this, it can be mentioned 
that the reflected populations play only a very small role in 
the main results of this paper. 
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Table 1. Comparison With Lemaire [1972] 



n p0 , % 

Lemaire [1972] 

This study 

0 + -H + crossover 

<0.0094 

5800 km 

4300 km 

ft,' at crossover 

£0.0094 

4.6 cm' 3 

25 cm' 3 

rc H + (10,000 km) 

<0.0094 

1.7 cm' 3 

5.9 cm' 3 

riQ+ (10,000 km) 

0.0 

8.0x1 O' 2 cm' 3 

2.7xl0' 2 cm' 3 

«o + (10,000 km) 

0.0094 

0.14 cm' 3 

2.9xl0" 2 cm* 3 

^ (10,000 km) 

0.0094 

2.2x1 O' 2 cm' 3 

3.6xl0' 2 cm' 3 

w H + (/j->oo) 

0.0094 

21 kins' 1 

21 km s' 1 


7.2. Comparisons With Previous Studies 

It is useful to compare our polar wind results with other cal- 
culations for this region. The first comparison will be with 
the results of Lemaire [1972], who calculated the steady state 
polar wind kinetically for all four species (0 + , H + , and thermal 
and superthermal electrons). This model was the first to con- 
sider photoelectrons in the polar wind, and predicted modest 
changes due to the new electron population. The low altitude 
boundary of Lemaire [1972] was at 950 km, with 7xl0 3 cm' 3 
for 0 + and 320 cm' 3 for H + , varying the photoelectron content 
from 0.0% to 0.0094%. These photoelectron concentrations 
are below those shown in section 7.1 but are valid composi- 
tions since the solar zenith angle is large at polar cap lati- 
tudes. The thermal plasma temperature was taken to be 3000 K 
everywhere, and the photoelectrons were assumed to have a 
characteristic energy of 10 eV. 

Results of this comparison are given in Table 1. As in the 
work of Lemaire [1972], our model also shows that n H + is not 
affected by the presence of photoelectrons at the concentra- 
tions taken for these calculations (0.0% and 0.0094%). Also, 


both models predict a crossover of the 0 + and H + densities at 
an altitude of less than an /?£. It is clear then that these con- 
centrations are very small and the polar wind could still be 
considered classical. Notice that in Table 1, the first three 
parameters shown have no significant dependence on the vary- 
ing photoelectron concentrations used by Lemaire [1972]. 
Also note that the models predict the same H + flow velocity at 
h -> «> of 21 km s' 1 . The other parameters, although quantita- 
tively different, show the same general trend, including a 
slight increase in the 0 + density at higher altitudes but not at 
the 0 + -H + crossover altitude. 

Since the recently developed model of Tam et ai [1995a] is 
also a steady state polar wind solution, this provides a good 
opportunity for assessment. It must be mentioned that the 
models are different since Tam et ai [1995a] is a numerical 
solution that includes collisions. We believe, however, that 
choosing the low-altitude boundary at 500 km, above the 
region of the source, is a more critical step than omitting col- 
lisional processes in the simulation region, becuase this pre- 
determines the solution above the lower boundary level [T. I. 
Gombosi, private communication, 1992]. Here we will point 
out the similarities and differences in the results. 

For this comparison, we tried to match the Tam et ai 
[1995a] boundary conditions as closely as possible. The 
plasma densities at 500 km were taken to be 6xl0 4 cm' 3 for 
0 + , IxlO 3 cm -3 for H + , and 0.1% of the balancing electrons for 
photoelectrons. The characteristic temperatures at the base are 
4500 K for the ions, 5000 K for the electrons, and 20 eV for 
the photoelectrons. Our simulation domain extended up to 5 
/?£, but only the results up to 2 R E are shown in the compari- 
son, as was presented by Tam et al. [1995a]. As pointed out in 
section 7.2, the location of the upper boundary has little effect 
on the results below 2 /?£. 

Figure 2 shows the results from the two models. In each 
plot, the photoelectron content used to produce each curve is 
indicated for our results, and the results of Tam et al. [1995a] 
are indicated with T95 and curves in bold. Similar quantities 
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Figure 2. Comparison of results with those of Tam et al. [1995], showing (a) ion density, (b) ion bulk 
velocity, (c) electrostatic potential, and (d) thermal electron temperature from 500 km to 2 /?£. Our results are 
indicated with the base photoelectron density used and theirs are indicated by T95 and bold linestyles. Low- 
altitude boundary conditions are 0 + density of 6xl0 4 cm' 3 , H + density of IxlO 3 cm' 3 , T i( f = 4500 K, T e(? =5000 
K, and E p( j= 20 eV. 
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are shown with the same linestyle, that is, <p( 0.1%) and 
<p{ T95) are both solid lines, while <p( 0.0%) is a dashed line. 

There is a good deal of agreement in the results of the 
models. In Figure 2a, it appears that above 2000 km the H + 
densities have the same scale height, and above 4000 km, the 
0 + densities are similar in magnitude. Notice that both models 
predict that 0 + is always the dominant ion. It is clear from 
Figure 2b that both models predict a much higher bulk 
velocity for H + than for 0 + . The tendencies for each species is 
also quite similar. Figure 2c shows that the potentials calcu- 
lated by the models above 2000 km are quite close. For this 
plot, <p( 2 R e ) is taken to be zero. Also, the thermal electron 
temperatures are comparable above 6000 km (Figure 2d). 

The difference in the results at low altitudes (s<2000 km) is 
clearly seen in the larger potential calculated in the Tam et at. 
[1995a] model. This is probably related to the very large and 
unrealistic electron temperatures that the Tam et al [1995a] 
model predicts for these altitudes, as compared with Akebono 
satellite data given by Abe et al. [1993]. The Tam et at. 
[1995a] model also calculates increased ion acceleration and 
decreased ion densities in this region. 

Our results reproduce the Akebono temperature profiles of 
Abe et al. [1993] very well in the low-altitude range. Figure 3 
shows this comparison. The data is from April 28 and May 
10, 1991, as the satellite passed from the dayside to the night- 
side at progessively higher altitudes. They presumed the data 
fluctuations were latitudinal or local time variations in the 
ionospheric conditions. Shown with the data are three curves 
at different levels of n p0 (0.02%, 0.03%, and 0.06%), which 
are reasonable numbers for the photoelectron content based on 
calculations from the Khazanov and Liemohn [1995] model for 
the conditions of an illuminated polar cap. This comparison 
shows that this data can be reproduced by this model for realis- 
tic estimates of the photoelectron content at the base. 

Also of interest in Figure 2 is the comparison of our results 
with and without photoelectrons. Figure 2a shows that the 
proton density is unaffected by the presence of photoelec- 
trons, being driven primarily by the changing magnetic field. 
This conclusion was first mentioned by Lemaire [1972], and 
was also illustrated earlier in this section. 0 + , however, 
shows a definite decrease in density without photoelectrons, 
and the 0 + /H + density ratio is unity just above 2 R E . From 
Figure 2c, the photoelectrons appear to increase (p 0 by nearly 
a volt, but (p(0.I%) and (p{ 0.0%) have comparable spatial dis- 
tribution trends. Comparing the thermal electron temperature 
profiles with and without photoelectrons, it is clear that there 
is an influence, The result without photoelectrons actually 
decreases with altitude just a bit, while the T e (0.\%) nearly 
doubles between 500 km and 2 R E . 

One other result plotted here is our model with a source term 
included in the thermal electron temperature equation. This 
source term represents energy deposition through Coulomb 
collisions with the photoelectrons. Although this term is in- 
consistent with the collisionless nature of the model, it is 
useful to illustrate the effects of including electron-electron 
collisions. This source term is expressed through the simple 
expression a(/t e [Khazanov et al. , 1992], where a 0 represents 
the energy input rate to a single thermal electron from colli- 
sions with photoelectrons. Here a value of 5xlO* 4 eV s" 1 is 
used, a typical heat input for a dayside electron at midlatitudes. 
This quantity is an overestimation of the heating expected in 
the polar region, since it was determined from a case where two 
conjugate ionospheres contribute to the photoelectron con- 
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Figure 3. Comparison of results for several values of n p0 
with data from Abe et al. [1993] from April 28 and May 10, 
1991. 


tent. This term replaces the zero on the right hand side of 
(26), and so contributes to the temperature equation (27) as 

r sUb 

J (a 0 n e /B)ds 

inside the integral on the right-hand side. The influence of 
this term is shown in Figure 2d, and it is clear that this source 
is not a significant factor in changing the thermal electron 
temperature in the polar wind, raising it by less than 1% 
throughout the altitude range. 

The issue of the calculation of the electron temperature is 
the most critical point of this paper. We used a fluid approach 
for the thermal electrons because our low-altitude boundary is 
at 500 km, where the thermal electron distribution function is 
very close to a Maxwellian and the thermal conductivity pro- 
cess plays a very important role in the formation of T e . 

In the work of Tam et al. [1995a], the polarization electric 
field expression allows for collisionless effects, and it 
depends on both T e \\ and T e ± [see, e.g., Demars and Schunk, 
1992, equation (28)]. In addition, Tam et al. [1995a] used the 
16-moment transport equations to describe the flow of both 
parallel and perpendicular thermal electron energy (through II 
and ± heat flows). In the present model, the Spitzer thermal 
conductivity was used to calculate the electron temperature, 
which requires a collision-dominated thermal electron popula- 
tion to be valid. However, the thermal electrons are not colli- 
sion-dominated in the polar wind above altitudes of about 
2000-3000 km [e.g., Schunk and Watkins, 1981; Demars and 
Schunk, 1989]. If the Spitzer conductivity is used for a colli- 
sionless plasma, erroneous T e profiles can result. Typically, 
the Spitzer conductivity tends to keep T e low by rapidly dis- 
tributing the heat along the geomagnetic field. This may 
partly explain the difference with the Tam et al. [1995a] calcu- 
lations, although we think that our treatment should not pro- 
duce a big difference with a more accurate treatment of the elec- 
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tron temperature. In fact, when considering only photoelec- 
tron heating of the thermal electrons, the difference in T e 
between noncollision-dominated and collision-dominated ap- 
proaches is only 25% [Gorbachev et ai, 1991]. Also, there is 
another comparison by Olsen and Leer [1996] for the solar 
wind between a 5-moment and an 8-moment description, and 
the results show a similar difference in the plasma temperature 
as was found by Gorbachev et al [1991]. 

8. Conclusion 

In order to answer questions about the maximum expected 
influence of photoelectrons on the polar wind, we developed a 
collisionless analytical model of the steady state polar wind 
with a kinetic description of the ions and photoelectrons and a 
fluid description of the thermal electrons. This model solves 
for the electrostatic potential and the electron temperature and 
describes the quantitative relation between the thermal plasma 
and the photoelectron parameters. The equations were itera- 
tively solved, using the conditions of quasi-neutrality and cur- 
rentlessness, and results of this study were presented. 

For low-altitude boundary temperatures of 2000 K for the 
ions and 2500 K for the electrons, varying the photoelectron 
content at the base produced significant changes in the elec- 
trostatic potential and thermal electron temperature. For 
n p( f= 0.01%, the potential drop is 2.87 V and the temperature 
reaches 3600 K; for n p Q=0.\0%, <po is 6.47 V with a high alti- 
tude temperature of 8800 K; and for n p o= 1.00%, (p 0 =l. 25 V 
and the temperature rises to 17,300 K. Increasing the low-alti- 
tude characteristic energies of the electron species increases 
the potential and electron temperature, while increasing the 
ion energies decreases the potential. Electron temperatures 
increase steadily with photoelectron content, always having 
the same tendency to rise quickly at low altitudes because of 
the electrostatic coupling and then reach an asymptotic tem- 
perature by 2 /?£. The results showed little influence on the 
location of the upper boundary, changing less than 5% in T e 
find less than 1% in <p when the boundary is moved from 5 R E 
to 3 R e . Inclusion of a source term to represent heating of the 
thermal electrons through Coulomb collisions with the pho- 
toelectrons increases the temperature by less than 1%. In the 
Ion results, a clear separation occurs between a classical polar 
wind at low photoelectron concentrations and an 0 + dominated 
polar wind for high photoelectron densities. 

The model was compared with results from several previous 
polar wind studies with photoelectrons. A comparison with 
that of Lemaire [1972] shows good agreement in the features 
of the results. Lemaire [1972], however, only considered pho- 
toelectron concentrations of less than a hundredth of a percent 
at the base, and thus the results are still within the classical 
polar wind description of supersonic H + dominance at high al- 
titudes. A comparison with the results of Tam et ai [1995a] 
shows that the two models show agreement at high altitudes 
but disagree near the lower boundary. We failed to reproduce 
the localized enhancement of electron temperature at 500- 
1500 km found by Tam et ai [1995a], but we have very good 
agreement with the electron temperature profile seen by the 
Akebono satellite at these altitudes [Abe et ai , 1993] for rea- 
sonable photoelectron concentrations as predicted by the 
model of Khazanov and Liemohn [1995]. 

The critical issue is the electron temperature calculation. 
This model uses a Spitzer conducitivity description to obtain 
T e , which theoretically is only valid for the thermal electrons 


in a collision-dominated region. This would only be true at 
altitudes less that 3000 km [Demars and Schunk , 1989]. 
However, it was shown by Gorbachev et al. [1991] and Olsen 
and Leer [1996] that this assumption should not be greatly dif- 
ferent from a more complicated description, and the model 
shows good agreement with the polar wind measurements of 
Abe et ai [1993], which indirectly demonstrates the validity 
of this assumption. 
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